cd  "E:\Replication Data for Media Freedom and Public Trust in Science" 

use Data_Gu_WGM2018, clear

**# 前期准备
***生成控制变量的宏
global cov ( i.female i.emp  i.urban_rural  i.income )
global agg ( c.loggdppc c.logpop c.log_nobel_pc c.log_qs500_pc c.e_pelifeex )

global cov2 (female emp urban_rural income)
global agg2 (loggdppc logpop log_nobel_pc log_qs500_pc e_pelifeex)

global control ///
female emp urban_rural income loggdppc logpop log_nobel_pc log_qs500_pc e_pelifeex

***画图设置
set scheme s2color
grstyle init
grstyle color background white
grstyle set plain,horizontal grid
grstyle set color Set1
grstyle color ci_area gs12%40

**#TABLE 1. Multilevel Model from the WGM 2018: Media Freedom Exposure
eststo clear

eststo m1: mixed y2 || country:, mle
estat icc
estadd scalar ICC = r(icc2)*100 

eststo m2: mixed y2 $cov2 edu $agg2 || country:, mle
estat icc
estadd scalar ICC = r(icc2)*100 

local var Mexp_conVDem_14after
eststo m3: mixed y2 $cov2 edu $agg2 `var' || country:, mle
estat icc
estadd scalar ICC = r(icc2)*100 
qui tab country if `var' != .
estadd scalar N2 = r(r)
eststo m4: mixed y2 $cov2 edu $agg2 c.`var'##c.edu || country:, mle
estat icc
estadd scalar ICC = r(icc2)*100 
qui tab country if `var' != .
estadd scalar N2 = r(r)

local var Mexp_binVDem_14after
eststo m5: mixed y2 $cov2 edu $agg2 `var' || country:, mle
estat icc
estadd scalar ICC = r(icc2)*100 
qui tab country if `var' != .
estadd scalar N2 = r(r)
eststo m6: mixed y2 $cov2 edu $agg2 c.`var'##c.edu || country:, mle
estat icc
estadd scalar ICC = r(icc2)*100 
qui tab country if `var' != .
estadd scalar N2 = r(r)

esttab m* using TG1.rtf , b(3) se(3) nonote nobaselevel noomitted  eqlabel(none) label  replace  star(+ 0.1 * 0.05 ** 0.01 *** 0.001) nomtitles  /// 
keep( _cons female emp urban_rural income edu loggdppc logpop log_nobel_pc log_qs500_pc e_pelifeex Mexp_conVDem_14after Mexp_binVDem_14after c.Mexp_conVDem_14after#c.edu c.Mexp_binVDem_14after#c.edu )  ///
order(_cons) ///
varlabel( _cons "Constant" female "Female" emp "Employment status" urban_rural "Residence" income "Income" edu "Education level" loggdppc "Log GDP per capita" logpop "Log population" log_nobel_pc "Log number of Nobel Prize winners per capita" log_qs500_pc "Log number of QS 500 universities per capita" e_pelifeex "Life expectancy"  Mexp_conVDem_14after "Continuous media freedom exposure" c.Mexp_conVDem_14after#c.edu "Continuous media freedom exposure × Education level" Mexp_binVDem_14after "Binary media freedom exposure" c.Mexp_binVDem_14after#c.edu "Binary media freedom exposure × Education level"  ) ///
mgroup("{\i DV = Trust in Science (IRT, from the WGM 2018)}") ///
stats( ICC N2 N, labels( "ICC (%)" "Country-Level Observations" "Individual-Level Observations") fmt( 2 0 0 ))  ///
title({\b Multilevel Model from the WGM 2018})  ///
note("{\i \b Notes:} + {\i p} < 0.1, * {\i p} < 0.05, ** {\i p} < 0.01,*** {\i p} < 0.001 (two-tailed test).")

**#TABLE 2: Fixed Effects Model from the WGM 2018: Media Freedom Exposure and Education Level
eststo clear

local var Mexp_conVDem_14after
reghdfe y2     ///
c.`var'##i.edu     ,cluster(country)  keepsing
local r2=e(r2_a)
eststo m1: expint y2  `var'  beta
estadd scalar r2a=`r2'
lincom primary-secondary
lincom primary-tertiary
lincom secondary-tertiary

local var Mexp_conVDem_14after
reghdfe y2     ///
c.`var'##i.edu     ,a(country ) cluster(country)  keepsing
local r2=e(r2_a)
eststo m2: expint y2  `var'  beta
estadd local rfe "√"
estadd scalar r2a=`r2'
lincom primary-secondary
lincom primary-tertiary
lincom secondary-tertiary

local var Mexp_conVDem_14after
reghdfe y2     ///
c.`var'##i.edu $agg##i.edu ,a(country i.birthyear) cluster(country)  keepsing
local r2=e(r2_a)
eststo m3: expint y2  `var'  beta
estadd local rfe "√"
estadd local bfe "√"
estadd local ctrl1 "√"
estadd scalar r2a=`r2'
lincom primary-secondary
lincom primary-tertiary
lincom secondary-tertiary

local var Mexp_conVDem_14after
qui reghdfe y2     ///
c.`var'##i.edu  $agg##i.edu $cov##c.`var'  ,a(country i.birthyear) cluster(country)  keepsing
local r2=e(r2_a)
eststo m4: expint y2  `var'  beta
estadd local ctrl2 "√"
estadd local ctrl1 "√"
estadd local rfe "√"
estadd local bfe "√"
estadd scalar r2a=`r2'
lincom primary-secondary
lincom primary-tertiary
lincom secondary-tertiary

local var Mexp_binVDem_14after
reghdfe y2     ///
c.`var'##i.edu $agg##i.edu $cov##c.`var' ,a(country i.birthyear) cluster(country)  keepsing
local r2=e(r2_a)
eststo m5: expint y2  `var'  beta
estadd local ctrl2 "√"
estadd local ctrl1 "√"
estadd local rfe "√"
estadd local bfe "√"
estadd scalar r2a=`r2'
lincom primary-secondary
lincom primary-tertiary
lincom secondary-tertiary

esttab m* using TG2.rtf , b(3) se(3) nonote nobaselevel  unstack noomitted  eqlabel(none) label  replace  star(+ 0.1 * 0.05 ** 0.01 *** 0.01) refcat(primary "{\i Effect of cumulative exposure to media freedom on}" , nolabel) /// 
varlabel(primary "Primary education or below" secondary "Secondary education" tertiary "College education or above") ///
mgroup("{\i DV = Trust in Science (IRT)}", pattern(1 0 0 0 0 0 0) ) ///
mtitles("{\i Continuous Media Freedom Exposure (V-Dem)}" "{\i Continuous Media Freedom Exposure (V-Dem)}" "{\i Continuous Media Freedom Exposure (V-Dem)}" "{\i Continuous Media Freedom Exposure (V-Dem)}" "{\i Binary Media Freedom Exposure (V-Dem)}") ///
stats(rfe bfe ctrl1 ctrl2 r2a N, labels("Country FE" "Birth year FE" "Individual-level controls × Media freedom" "Country-level controls × Education level" "Adjusted R²"  "Observations") fmt( 0 0 0 0 3 0 ) ) ///
title({\b Fixed Effects Model from the WGM 2018: Media Freedom Exposure and Education Level}) ///
note("{\i \b Notes:} + {\i p} < 0.1, * {\i p} < 0.05, ** {\i p} < 0.01,*** {\i p} < 0.001 (two-tailed test)." )








